

#include "CurrentPath.h"

namespace falko
{



Loop::Loop(hpv3d_t const &center, hpf_t const &radius) : BaseCurrentPath()
{
	path_center_ = center;
	radius_ = radius;
	partitions_ = 50;
}

void Loop::getField(hpv3d_t const &point, hpv3d_t &field)
{
	if (inContact(point))
		return;
	PathFactory::solve(&field, &point, this);
}

bool Loop::inContact(hpv3d_t const &point)
{
	return false;
}

void Loop::getPartitionLength(hpf_t &Length) const
{
	Length = ((_PI * 2) * radius_) / partitions_;
}

void Loop::getPoint(uint const i, hpv3d_t &point) const
{
//		HpFileWriter sout;
	hpf_t theta = ((_PI * 2) / partitions_ ) * i;
//		sout << "theta" << "..........." << theta << "\n";

	hpf_t cos_theta;
	hpf_t sin_theta;

	hpf_cos(cos_theta, theta);
//		sout << "cos_theta" << "..........." << cos_theta << "\n";

	hpf_sin(sin_theta, theta);
//		sout << "sin_theta" << "..........." << sin_theta << "\n";

	point.x_ = path_center_.x_ + radius_ * cos_theta;
//		sout << "point.x" << "..........." << point.x_ << "\n";

	point.y_ = path_center_.y_ + radius_ * sin_theta;
//		sout << "point.y" << "..........." << point.y_ << "\n";

	point.z_ = path_center_.z_;
//		sout << "point.z_" << "..........." << point.z_ << "\n";
}


MagneticField::~MagneticField()
{
	for(int i = paths_.size()-1; i >= 0; --i)
	{
		if (paths_[i].own_)
			delete paths_[i].path_;
	}
}

void MagneticField::addPath(BaseCurrentPath* path, bool own)
{
	paths_.push_back(PathPointer(path, own));
}

void MagneticField::getField(hpv3d_t const &point, hpv3d_t &field) const
{
	field = 0.0;
	hpv3d_t temp;
	for(int i = paths_.size()-1; i >= 0; --i)
	{
//		temp = 0.0;
		paths_[i].path_->getField(point, temp);
		field += temp;
	}
}

bool MagneticField::inContact(hpv3d_t const &point)
{
	return false;
}





void PathFactory::MakeCoil(MagneticField *path,
						   uint const loops,
						   hpf_t const &wire_radius,
						   hpf_t const &length,
						   hpf_t const &coil_inner_radius,
						   hpf_t const &current,
						   hpv3d_t const &center)
{
	uint n_length = static_cast<uint>(ceil(((length/(wire_radius * 2)).toDouble())));
	hpv3d_t start_center = center;
	start_center.z_ += length/2;
	hpv3d_t init_center = start_center;
	hpf_t init_radius = coil_inner_radius;

	for (uint i = 0; i < loops; ++i)
	{
		path->addPath(new Loop(init_center, init_radius, current), true);

		if (i % n_length)
		{
			start_center.z_ += wire_radius * 2;
		}
		else
		{
			init_center = start_center;
			init_radius += wire_radius * 2;
		}
	}

}


void PathFactory::solve(hpv3d_t *field, hpv3d_t const *point, BaseCurrentPath const *path)
{
	hpf_t   dlLength(0.0);
	hpv3d_t dl(0.0, 0.0, 0.0);
	hpv3d_t point_2(0.0, 0.0, 0.0);
	hpv3d_t displacement(0.0, 0.0, 0.0);
	hpv3d_t intergrand(0.0, 0.0, 0.0);
	*field = 0.0;
//	HpFileWriter sout;
	uint partitions = path->getPartitions();
	for (uint i = 0; i < partitions ; ++i)
	{

		path->getPoint(i, dl);
//		sout << "path->getPoint(i, dl)" << "..........." << dl << "\n";


		path->getPoint(i+1, point_2);
//		sout << "path->getPoint(i+1, point_2);" << "..........." << point_2 << "\n";

		dl = point_2 - dl;
//		sout << "dl = point_2 - dl;" << "..........." << dl << "\n";

		dl.normalize();
//		sout << "dl.normalize();" << "..........." << dl << "\n";

		displacement = *point - point_2;
//		sout << "displacement = *point - point_2;" << "..........." << displacement << "\n";

		intergrand =  dl % displacement;
//		sout << "intergrand =  dl % displacement;" << "..........." << intergrand << "\n";

		intergrand /= displacement.magnitude().pow(3);
//		sout << "intergrand /= displacement.magnitude().pow(3);" << "..........." << intergrand << "\n";

		*field += intergrand;
//		sout << "*field += intergrand;" << "..........." << *field << "\n";
	}
	path->getPartitionLength(dlLength);
	*field *= dlLength;
	*field *= (PERM_FREE_SPACE * (path->current())) / (_PI * 4.0);
}


}
